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Abstract. Star clusters with a high central density con¬ 
tain an ecological network of evolving binaries, affected 
by interactions with passing stars, while in turn affecting 
the energy budget of the cluster as a whole by giving off 
binding energy. This is the first paper in a series aimed at 
providing the tools for increasingly realistic simulations of 
these ecological networks. 

Here we model the core of a globular star cluster. The 
two main approximations are: a density of stars constant 
in space and time, and a purely single star population in 
which collisions between the evolving stars are modeled. 
In future papers in this series, we will relax these crude 
approximations. Here, however, they serve to set the stage 
before proceeding to the additional complexity of binary 
star interactions, in paper H, and background dynamical 
evolution, in later papers. 


Key words: methods: numeric - celestial mechanics: stel¬ 
lar dynamics - stars: evolution - stars: blue stragglers - 
globular clusters: general 

1. Introduction 

In dense stellar systems, such as open and globular clusters 
and galactic nuclei, encounters between individual stars 
and binaries can affect the dynamical evolution of the sys¬ 
tem as a whole on a time scale comparable to, or even 
shorter than, a Hubble time. In order to reach a detailed 
theoretical understanding of such systems, the following 
three steps are necessary. 

First, we need to understand the basic mechanism of 
the dynamical evolution, in the limit of a point-mass ap¬ 
proximation for the stars. Second, effects of dynamical en¬ 
counters on the internal evolution of single stars and bina- 
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ries has to be taken into account. Third, we have to model 
the feedback of these internal changes onto the dynamical 
evolution of the whole system. Let us briefly review each 
step. 

Great progress has been made with the first step, mod¬ 
eling the dynamical evolution of point-mass systems. In 
the seventies, the processes of core collapse and mass 
segregation were studied with the use of various types 
of Fokker-Planck approximations. In the eighties, these 
simulations were extended successfully beyond core col¬ 
lapse, and various studies were made of the phenomenon of 
gravothermal oscillations, ubiquitous in the post-collapse 
phase. Of these models a few even include mass loss due to 
the evolution of the stars (Chernoff & Weinberg 1990). In 
the nineties, we are finally beginning to switch over from 
Fokker-Planck approximations to much more detailed and 
realistic iV-body simulations. In 1995, the construction of 
the GRAPE-4, a special-purpose machine with Teraflops 
speed, has made a 32, 000-body simulation feasible, pro¬ 
viding the first direct evidence of gravothermal oscillations 
(Makino, 1996a,b). Extending these simulations to the full 
realm of globular clusters {N = 10^ ~ 10®) will require 
Petaflops speed, something that could be realized by fu¬ 
ture special-purpose machines in the GRAPE series by as 
early as the year 2000. 

While the point-mass approximation provides a good 
qualitative guide for the construction of dynamical models 
of dense stellar systems, this approximation quickly breaks 
down when we require quantitatively accurate results. The 
second step attempts to model the effects of close encoun¬ 
ters. A number of different investigations have estimated 
the rate at which physical collisions have take place, un¬ 
der various circumstances (Hills & Day 1976, Verbunt & 
Meylan 1988, Di Stefano & Rappaport 1992 and Davies 
& Benz 1995). However, little progress has been made so 
far in following the changes induced in the stellar popula¬ 
tion, beyond enumerating the number of mergers. In the 
simulations presented below, collisions are modeled in an 
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evolving population of single stars in a high-density stellar 
environment. Paper II in this series will extent our treat¬ 
ment to follow the induced changes in binary systems, 
both on the level of changes in orbital parameters as well 
as in the internal structure of the stars. 

However, such investigations are only a start, and can¬ 
not lead to a quantitative modeling of dense stellar sys¬ 
tems, since they are not yet self-consistent. What is needed 
in addition is a treatment of the feedback mechanism, from 
the changes in single stars and binaries back to the overall 
dynamics of the system. This third step is being pioneered 
for open clusters by Aarseth (1996). The current series of 
papers aims to provide self-consistent models of this type, 
by coupling relatively crude stellar evolution recipes, doc¬ 
umented and tested in papers I and II, to a fully dynamical 
Wbody system. 

This paper is organized as follows. Our approach to 
the study of the ecology of star clusters is summarized 
in somewhat more detail in §2. The next section, §3, de¬ 
scribes our simulation techniques, and the various approx¬ 
imations involved. In §4, we present the result of a sim¬ 
ulation starting a single model, with a minimum of free 
parameters. The results of a more realistic core model run 
are presented in §5. §6 sums up. 

2. Star Cluster Ecology 

Stellar evolution plays a role in star cluster evolution sim¬ 
ilar to the role played by nuclear physics in stellar evo¬ 
lution. In both cases, the microphysical processes play a 
crucial role in the mechanism of energy generation in the 
central parts of the system under consideration, a mecha¬ 
nism that tries to balance the energy losses at the outskirts 
(tidal radius and photosphere, respectively). 

In the next two subsections we look separately at the 
different forms of physics input necessary to follow the 
evolution of a star cluster. The third subsection then dis¬ 
cusses their interconnection. For future reference, initial 
conditions are discussed in the fourth subsection, while 
the last subsection provides a brief outline of the series of 
papers of which this one is the first. 

2.1. Stellar Dynamics Simulations 

Great progress has been made in the study of star cluster 
dynamics, using various approximate methods in which 
the stars have been treated like a form of fluid, either 
three-dimensional as in conducting gas sphere models, or 
six-dimensional as in Fokker-Planck models. In both cases, 
the main effect of encounters has been taken into account 
by a form of effective two-body relaxation. We refer to 
Hut et al. (1992) for a review of these methods. 

Unfortunately, both methods have two intrinsic hand¬ 
icaps that make them unsuitable for a detailed quantita¬ 
tive modeling of the evolution of a globular cluster past 
core collapse. First, they are not set up to deal with the 


separate evolution of internal and external degrees of free¬ 
dom of the binaries that play a central role in the energy 
generation processes in the cluster. 

The second problem stems from an introduction of 
a mass spectrum, as well as a distinction between stars 
of different radii, such as dwarfs, main-sequence stars, 
and giants. The root of the problem here is that a gas 
sphere or Fokker-Planck approach does not follow indi¬ 
vidual stars, but rather distribution functions. When the 
number of independent parameters characterizing the dis¬ 
tribution functions becomes too large, there will be less 
than one star left in a typical cell in parameter space — 
something that clearly invalidates the statistical hypoth¬ 
esis on which these methods are based. 

The only solution is to drop the statistical assumption, 
and to revert to a star-by-star modeling of a globular clus¬ 
ter, through direct iV-body calculations. The draw back 
of such an approach has long been the prohibitive calcu- 
lational costs involved, and until recently typical produc¬ 
tion runs only included a few thousand stars. To extend 
such numbers to include several hundred thousand stars, 
characteristic of realistic globular clusters, requires an in¬ 
crease of two orders of magnitude in star number, or a 
factor million in computational cost, from Gigaflops days 
to Petaflops days (Hut et al. 1988). 

Recently, the number of stars modeled in direct N- 
body calculations has been increased significantly, to iV = 
32, 000, using the GRAPE-4, a form of special-purpose 
hardware developed by a group of astrophysicists at Tokyo 
University, running at a speed of I Tflops (Makino 1996a). 
The first scientific results of the GRAPE-4, including the 
first convincing evidence of gravothermal oscillations in 
fV-body simulations, predicted by Sugimoto & Bettwieser 
(1983), have been presented by Makino (1996ab). 

The next, and definitive step that will enable any glob¬ 
ular cluster to be modeled realistically might take place 
as early as the year 2000. If funding can be found, there is 
no technological obstacle standing in the way of a speedup 
of the current GRAPE-4 machine by a factor of a thou¬ 
sand, during the next five years. Most of this speed-up 
will come from further miniaturization, allowing a larger 
number of gates to be mounted on a single chip, and al¬ 
lowing a higher clock speed as well. A Petaflops machine 
by the year 2000, allowing simulations of core collapse and 
post-collapse evolution with up to 10® particles, is thus a 
realistic goal. 

2.2. Stellar Evolution Population Synthesis 

The first serious attempts to understand and simulate the 
evolution of close binaries were made in the mid thirties 
(Haffner & Heckmann 1937) and late fifties by Grawford 
(1955), Kopal (1956) and Huang (1956) followed by Mor¬ 
ton (1960) and the standard work in binary evolution from 
Kippenhahn & Weigert (1967). Synthesis of complete pop¬ 
ulations of single stars became popular in the mid seven- 
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ties when Tinsley & Gunn (1976) simulated the giant- 
branch luminosity functions for giant elliptical galaxies. 
However, it is only recently that detailed studies simu¬ 
late complete populations of close binaries starting with 
Dewey & Cordes (1987) who tried to understand the evo¬ 
lutionary sequence of radio pulsars and the presence of 
an asymmetry in the velocity distribution of single ra¬ 
dio pulsars. In later papers, similar evolutionary scenarios 
for the formation of binary neutron stars were studied in 
more detail (Tutukov & Yungelson 1993, Lipunov et al. 
1995, Portegies Zwart & Spreeuw 1996 and Lipunov et al. 
1996), for high mass X-ray binaries and the supernova 
rate in the galaxy (Tutukov et al. 1992, Lipunov 1994, 
Dalton & Sarazin 1995, Portegies Zwart & Verbunt 1996) 
and for lower mass systems with a neutron star (Web- 
bink & Kalogera 1993, Pols & Marinus 1994) or a white 
dwarf (de Kool 1990, Kolb & Ritter 1992) as the accret¬ 
ing object. The Hrst couragous attempt to combine stellar 
and binary evolution within the collisional evironment of 
a globular cluster was performed by Sutantyo (1975), fol¬ 
lowed more recently by Di Stefano & Rappaport (1992), 
Sigurdsson & Phinney (1993), Leonard (1994), Davies & 
Benz (1995) and Davies (1995). 

2.3. Ecological Networks 

Purely stellar-dynamical calculations often rely on rather 
severe approximations, such as a representation of stars 
by equal-mass point-masses. And there is a good reason 
for doing so, since any single deviation from that simple 
recipe requires other deviations as well. Let us look at one 
example. 

As soon as we introduce a mass spectrum in a star 
cluster simulation, we will see that the heavier stars start 
sinking toward the center, on the dynamical friction time 
scale, shorter than the two-body relaxation time by a fac¬ 
tor proportional to the mass ratio of individual heavy stars 
with respect to that of typical stellar masses. The reason 
is that relaxation tends toward equipartition of energy, 
which implies that heavier stars will move more slowly 
and therefore gather at the bottom of the cluster poten¬ 
tial well. 

If stars would live forever, there would be a large over¬ 
concentration of heavy stars in the core of a star cluster. 
However, in reality there is an important counter-effect: 
heavy stars burn up much faster than lighter ones. They 
may or may not leave degenerate remnants, that may 
or may not be heavier than the average stellar mass in 
the cluster (a quantity that decreases in time). Clearly, it 
would be grossly unrealistic to introduce a mass spectrum 
without removing most of the mass of the heaviest stars 
on the time scale of their evolution off the main sequence 
and across the giant branches. 

Another reason for introducing finite life times for stars 
comes from abandoning the very restrictive point mass 
model. As soon as we do that, giving our stars a Hnite 


radius will give rise to stellar collisions. The heavier stars 
produced in the collision of two turn-off stars, for exam¬ 
ple, will burn up on a time scale an order of magnitude 
smaller than the age of the cluster. Again, we have to 
take this into account to be consistent, especially since 
the merger products themselves are prime candidates for 
further merging collisions. 

The need to let many stars shed most of their mass, to¬ 
gether with the fact that most of the energy in a globular 
cluster is locked up in binaries, poses a formidable consis¬ 
tency problem. Since binaries play a central role in cluster 
dynamics, consistency requires that we follow their com¬ 
plex stellar evolution, which involves mass overflow (which 
can be stable or unstable, and can take place on dynamical 
or thermal or nuclear time scales) and the possibility of a 
phase of common-envelope evolution. On top of all that, 
we will have to Hnd simple recipes for the hydrodynamic 
effects occurring in three-body and four-body reactions, 
and in occasional IV > 4 reactions, which are bound to 
occur in dense cluster centers. 

To sum up: there does not seem to be a half-way stop¬ 
ping point, at which we can expect to carry out consistent 
cluster evolution simulations. Either we study the inter¬ 
esting but unrealistic mathematical-physics problem of an 
equal-mass point particle model, or we opt for a realistic 
model with some set of stellar-evolution recipes. The main 
question here is: what is the simplest set that is still con¬ 
sistent? 

2 . 4 . Initial Conditions 

In most stellar dynamics simulations of star clusters, the 
Plummer model is used as a standard model to specify the 
initial conditions for the distribution of the point particles. 
While not very realistic, this choice has had the advantage 
of making comparisons between different runs, as well as 
between different approaches, relatively straightforward. 
Of course, when attempts are made to model particular 
star clusters, other models have often replaced the Plum¬ 
mer model as a starting point. King models, for example, 
are already more realistic in that they provide a form of 
spatial cut-off that can be interpreted as a tidal radius. 

For similar reasons, we will use a standard model for 
our simulations that combine stellar dynamics and stellar 
evolution. In most cases, the use of our standard model 
will be mostly for illustrative reasons, to provide a gauge 
for comparison between our various results, as well as be¬ 
tween our results and that of others. For historical reasons, 
we choose our standard model to be based on a Plummer 
model for the macroscopic initial star distribution, and a 
Salpeter model for the initial mass function. 

An additional advantage of these simple choices is that 
they limit the number of free parameters. The Plummer 
model, for example, contains only one free parameter, N, 
the number of stars in the system (apart from a choice of 
mass and length scales, that are irrelevant in the point par- 
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tide case). In contrast to the Plummer model, our stan¬ 
dard model can be expected to form a multi-parameter 
family. As soon as we abandon the point-mass approach, 
we have to deal with microscopic as well as macroscopic 
mass and length scales. 

Of these various scales, the macroscopic quantities can 
be chosen independently, while the microscopic ones can 
be fixed, statistically, by specifying a mass distribution 
together with appropriate cut-off masses at the high and 
low end. In general, an arbitrary functional form for the 
mass distribution function can lead to an arbitrarily large 
number of parameters. Interestingly, our standard model 
definition allows us to limit the total number of free pa¬ 
rameters to three. 

Starting from the macroscopic side, we can take the 
total mass M and the half-mass radius Vh of the Plummer 
model as our first two free parameters. With a Salpeter 
choice of powerlaw distribution function, the third free 
parameter can be chosen in the form of the lower mass 
cut-off m_. The higher-mass cut-off to+ could be speci¬ 
fied independently, but this is not strictly necessary: since 
the Salpeter distribution function converges at the high- 
mass end, we can simply parcel out the total mass M over 
different stellar masses, between star masses of m = m_ 
and TO+ = oo, and we will naturally be left with a sin¬ 
gle most massive star. This procedure is not unrealistic: 
nature probably limits the number of high-mass stars in 
medium-size galactic clusters in a similar way. 

In fact, we can go even further, and make the follow¬ 
ing somewhat arbitrary but natural choices: = 10 pc, 

m_ = O.IM©. This leaves only the total mass M to be 
specified, or equivalently, the total number of stars N. 
For future convenience, we will refer to this ‘most stan¬ 
dard’ model as our reference model. For systems with a 
few thousand stars, we are dealing with a typical open 
cluster, with velocity dispersions of order 1 km/s, while for 
a few hundred thousand stars, we have a reasonable ap¬ 
proximation to a globular cluster, for which typical stellar 
velocities are an order of magnitude higher. 

In addition to this standard model, the various papers 
in this series will also contain the results of more realistic 
models. However, we will typically provide at least one run 
from a standard model, in order to provide comparison 
material for the more detailed models. 

2.5. Stepping Stones 

In the current series of papers, our goal is to provide a 
series of ecological simulations, based on a flexible stellar 
dynamics code coupled to a comprehensive set of stellar 
evolution modules. These modules in turn are based on 
recipes that govern the behavior of both single star and bi¬ 
nary star evolution, as well as interactions between larger 
numbers of stars. 

In order to present results that can be reproduced and 
critically assessed by other groups, we clearly document 


the recipes used, as well as their coupling to the dynam¬ 
ics. With this aim, we give a detailed description of our 
approach in the first few papers in this series, which will 
form stepping stones towards a full-fledged ecological star 
cluster evolution code. 

The present paper starts off with rather extreme ap¬ 
proximations for the stellar dynamics, as well as the stellar 
evolution parts of our simulations. With respect to the for¬ 
mer, we start with a laboratory-type situation, in which 
we consider a homogeneous distribution of stars, kept con¬ 
stant in time. With respect to the latter, we consider a 
population of single stars only. Paper II will relax the sec¬ 
ond assumption, by introducing a population of primor¬ 
dial binaries, and allowing the formation of new binaries 
as well. Later papers will subsequently relax the former as¬ 
sumption, with the ultimate goal of using a self-consistent 
Wbody code. 

3. A Static Homogeneous Environment with Sin¬ 
gle Stars 

3.1. Initial Conditions 

In the present paper, we keep the dynamical environment 
as simple as possible, in order to focus on the stellar evo¬ 
lution recipes, that are introduced here and used in sub¬ 
sequent papers as well. The stellar distribution is take to 
be in thermal equilibrium, with a density that is constant 
in space and time. In addition, an additional simplifica¬ 
tion is obtained by excluding any primordial binaries, and 
ignoring binary formation channels. Within this setting, 
random encounters between single stars will lead to col¬ 
lisions resulting in the formation of merger products, the 
evolution of which can then be followed along with the 
evolution of the original single stars. 

3.1.1. Initial Mass Function 

While our main aim is to set-up and clarify our stellar 
evolution recipes, we present two calculations that could 
be interpreted as having a limited astrophysical interpre¬ 
tation, one for the core of an w-Centauri-like cluster (§4), 
and one for the core of an M-15-like cluster (§5). Our 
choice of constant density implies that we can only hope 
to model the history of a cluster core, not that of a cluster 
as a whole. To specify the mass distribution, we first take 
our standard choice: a Salpeter initial mass function (§4), 
which we will use to model a relatively unevolved core. 
Our second choice will be a much more flat distribution, 
which is more appropriate for a high-density post-collapse 
cluster core (§5). 

3.1.2. Mass and Number Densities 

If we specify the mass density p for the stars in our clus¬ 
ter case, we can use the mass function to determine the 
number density n = n(p). In a homogeneous medium the 
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relation is linear, and for the simplest case of a powerlaw 
mass function /(m) oc m~°‘, we find 

n _ a-2 m]r°‘ - 
p a — 1 ’ 


where m_ and m+ are the lower and upper mass cut-offs, 
respectively. For the example of an initial Salpeter mass 
function, a = 2.35, we find 


n __ 0.26 
p m- ’ 


( 2 ) 


when we neglect the fact that the upper mass cut-off is fi¬ 
nite. The inverse quantity niav is the average stellar mass: 


TTiav = 3.9m_ = O.39M0 


( 3 ) 


in our standard case where we take a lower cut-off mass of 
O.IM©. The median mass rrimed for a Salpeter distribution 
is 


which we can then compare with that of an astrophysical 
system, such as the core of a globular cluster. 

For an equal-mass cluster model that is close to ther¬ 
modynamic equilibrium, the density drops by roughly a 
factor three, from the center to the edge of the core. This 
implies that the local density of collisions, which is pro¬ 
portional to the square of the density, drops by an order of 
magnitude. In the more realistic case of a mass spectrum 
the situation is even worse, since the density of the heav¬ 
ier stars drops off faster than that of the lightest stars. 
In the present paper we will not attempt to model these 
density dependent effects, and instead we will keep the 
density of all mass groups constant throughout the region 
of our simulation. It is clear, therefore, that our results are 
mainly for the purpose of illustration, and that any com¬ 
parison with actual systems will have to be taken with 
many grains of salt. 

The only question remaining is the definition of a core 
radius re- For an equal mass system, we have (Spitzer 
1987) 


TOmed = = O.17M0 (4) 

which means that in our standard population most stars 
have a mass well below O. 2 M 0 . 

Even in our simple case of a homogeneous system, the 
linear relationship n{p)/p = Ifmav involves a complicated 
time dependent factor mav Not only does the upper mass 
cut-off TO+ (roughly the main sequence turn-off mass) de¬ 
pend on time, but what is worse, the distribution of rem¬ 
nants, in the form of black holes, neutron stars and white 
dwarfs, does not obey any simple power law, even if their 
progenitors did. In general, therefore, the coefficient n/p 
has to be determined numerically, as a function of time. 

3.1.3. Velocity Dispersions 

In thermal equilibrium, equipartition of kinetic energy 
tells us how the velocity dispersions scale for stars with 
different masses. We only have to specify the three- 
dimensional velocity dispersion v for one particular mass, 
say u(lM 0 ) = V 0 , in order to determine the 3D velocity 
dispersion z;(m) for stars of general mass m: 

v{m) = - V 0 . (5) 

\ m J 

3.1.4. Core Radius and Core Mass 

The three choices discussed so far, namely that of an initial 
mass function, a density, and a temperature, specify the 
intensive thermodynamic properties. This in turn enables 
us to calculate the local rate of collisions, per unit time, 
and per unit volume. In order to extract global informa¬ 
tion, we have to specify extensive quantities as well, such 
as the total volume or total mass of our system. This will 
allow us to determine a global collision rate per unit time. 



( 6 ) 


In the presence of a mass spectrum, we have to modify this 
equation. Although the velocity dispersion is now quite 
different for different mass groups, the average kinetic en¬ 
ergy per star {l/N)Ekin is independent of mass, with N 
the total number of stars in the core. Rewriting the above 
formula, we have: 


3 2 

AiiGp M 


^kin 1 


( 7 ) 


where M is the core mass. For a general mass spectrum, 
we can substitute Ekin = (-N/ 2 )M 0 V 0 ^ which gives 


Tc = 


— ) 

47rG/9 ) 


1/2 


mav ) 



( 8 ) 


In this expression, the right hand side contains only local 
quantities, and the global quantity Tc is given in terms of 
those. 

Note that this is not the only possible generalization of 
the equal-mass expression, but it is a natural one, and it 
reverts to the original expression in cases where the mass 
of the core is dominated by stars in a relatively small mass 
range, as is the case, for example, in a post-collapse core 
of a globular cluster. 

Other global quantities can be derived from Vc, such 
as the core mass M: 

M = ^irrlp, (9) 


where we have used the fact that in an isothermal sphere 
the average density in the core is roughly half the central 
density, a relationship, while not exact, is certainly good 
enough for our purpose of relating our results to astro- 
physical systems. 
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3.2. Recipes for Stellar Evolution 

The stars in our computations are evolving. To describe 
their evolution, we use the formulae htted to the results 
of full stellar evolution calculations, by Eggleton et al. 
(1989). These formulae give the radius and luminosity 
(for population I stars) as a function of time, on the 
main sequence, in the Hertzsprung gap, the (sub)giant 
branch, on the horizontal branch, and on the asymptotic 
giant-branch. We use these population I recipes, because 
the more appropriate data for population II stars are not 
available in the same convenient form. In addition to the 
radius, we need the core mass for stars that have left the 
main sequence. We derive these from the luminosity, ac¬ 
cording to Eggleton et al. (1989), and core-mass luminos¬ 
ity relations, according to Boothroyd & Sackmann (1988), 
Paczyhski (1970) and Iben & Truran (1978). The details 
of this procedure are described in Portegies Zwart and 
Verbunt (1996, Section 2.1). 

3.3. Recipes for Individual Encounters 

In this paper, we only treat single stars, and accordingly 
the only outcome allowed for a close encounter is a merged 
object. The merging between the two stars in an encounter 
in our calculation is generally assumed to conserve mass, 
which in fact may be a reasonable approximation (Benz 
and Hills 1987, 1989, 1992, Rasio and Shapiro 1991, 1992). 

Only a limited number of simulations of encounters 
between stars has been performed, and these does not 
cover all possible combinations that may occur in a clus¬ 
ter. Also, different authors do not agree on the details of 
the outcomes for the same type of encounter. We there¬ 
fore have chosen to use a set of simple prescriptions for 
the outcome of stellar collisions, often chosen without de¬ 
tailed justification. In the future these prescriptions can 
be refined, when more accurate calculations for collisions 
become available. Meanwhile, our results will help in de¬ 
termining which of all possible types of encounter are most 
frequent, and therefore deserve closer attention. 

We describe our treatment of the possible outcomes of 
the encounters of two stars ordered by the evolutionary 
state of the more massive of the two, the primary. Table 1 
summarizes this treatment. 

3.3.1. Main-sequence primary 

If both stars involved in the encounter are main-sequence 
stars the less massive star is accreted conservatively onto 
the most massive star. The resulting star is a rejuvenated 
main-sequence star (see Lai et al. 1993, Lombardi et al. 
1995). The details of this procedure are described in Ap¬ 
pendix C4 of Portegies Zwart and Verbunt (1996). 

If the less massive star in the encounter has a well 
developed core (giant or subgiant) this core is treated as 
the core of the merger product. The main-sequence star 
and the envelope of the giant are added together to form 


Table 1. Simplified representation of possible merger out¬ 
comes. The four columns correspond to the four choices given 
for the type of massive star (primary), while the four rows in¬ 
dicate the type of less massive star (secondary): main-sequence 
star (ms), (sub)giant (sg), white dwarf (wd) and neutron star 
(ns). In this table we do not discriminate between stars in the 
Hertzsprung gap (Hg) or on the first and second ascent on the 
asymptotic-giant branch (AGB). 


star 

primary 

ms sg wd ns 

ms 

wd ns 

ms sg -1- -1- 

disc disc 

sg 

wd ns 

Hg AGB -b -b 

disc disc 

wd 

sg AGB 

ns 

TZO TZO 


the new envelope of the merger. In general the mass of 
the core is relatively small compared to its envelope and 
the star is assumed to continue its evolution through the 
Hertzsprung gap. Note that this type of encounter can 
only occur when the main-sequence star is in itself a col¬ 
lision product (e.g. a blue straggler). 

When a main-sequence star encounters a less massive 
white dwarf, we assume that the merger product is a gi¬ 
ant, whose core and envelope have the masses of the white 
dwarf and the main-sequence star, respectively. We then 
determine the evolutionary state of the merger product, 
as follows. We calculate the total time tagb that a sin¬ 
gle, unperturbed star with a mass equal to that of the 
merged star spends on the asymptotic giant-branch, and 
the mass Wc.agb of its core at the tip of the giant branch. 
The age of the merger product is then calculated by adding 
tagbWc/mc,agb to the age of an unperturbed star with 
the same mass at the bottom of the asymptotic giant 
branch. For example, a single, unperturbed 1.4M0 star 
leaves the main-sequence after 2.52 Gyr, spends 60Myr 
in the Hertzsprung gap, moves to the horizontal branch 
at 2.96 Gyr, and reaches the tip of the asymptotic gi¬ 
ant branch after 3.06 Gyr, with a core of 0.64 Mq. Thus, 
if a 0.6 Mq white dwarf mergers with an 0.8 Mq main- 
sequence star, the merger product has an age of 2.87 Gyr, 
leaving it another 180 Myr before it reaches the tip of the 
asymptotic giant-branch. 

If the less massive star is a neutron star a Thorne 
Zytkow object (Thorne & Zytkow 1977) is formed. 
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3.3.2. Evolved primary 


Collisions between white dwarfs are ignored. 


When a (sub)giant or asymptotic branch giant encoun¬ 
ters a less massive main-sequence star, the main-sequence 
star is added to the envelope of the giant, which stays 
in the same evolutionary state, i.e. remains a (sub)giant, 
c.q. asymptotic branch giant. Its age within that state 
is changed, however, according to the rejuvenation cal¬ 
culation described in Section C.3 of Portegies Zwart and 
Verbunt (1996). For example, an encounter of a giant of 
0.95 Mq and age 11.34 Gyr with a 0.45 Mq main-sequence 
star produces a giant of 1.4 Mq with an age of 2.67 Gyr. 

When both stars are (sub)giants the two cores are 
added together and form the core of the merger product 
(see also the results of the smoothed particle hydrody¬ 
namics computations performed by Davies et al. 1991 and 
Rasio and Shapiro 1995). Half the envelope mass of the 
(less massive) encountering star is accreted onto the pri¬ 
mary. The merger product continues its evolution starting 
at the next evolutionary state; thus a (sub)giant contin¬ 
ues its evolution on the horizontal branch and a horizontal 
branch star becomes a asymptotic-giant branch star. The 
reasoning behind this assumption is that an increased core 
mass corresponds to a later evolutionary stage. 

If the less massive star is a white dwarf then its mass 
is simply added to the core mass of the giant, and the 
envelope is retained. If the age of the giant before the en¬ 
counter exceeds the total life time of a single unperturbed 
star with the mass of the merger, then the newly formed 
giant immediately sheds its envelope, and its core turns 
into a single white dwarf; if not then the merged giant 
is assumed to have the same age (in years) as the giant 
before the collision, and continues its evolution as a single 
unperturbed star. 

If the encountering star is a less massive neutron-star 
a Thorne Zytkow object is formed. 


3.3.3. White-dwarf primary 

In an encounter between a white dwarf and a less massive 
main-sequence star, the latter is completely disrupted and 
forms a disk around the white dwarf (Ruffert & Muller 
1990, Rasio & Shapiro 1991). The white dwarf accretes 
from this disk at a rate of one percent of the Eddington 
limit. If the mass in the disc exceeds 5% of the mass of 
the white dwarf, the excess mass is expelled from the disc 
at a rate equal to the Eddington limit. 

If a white dwarf encounters a less massive (sub)giant, 
a new white dwarf is formed with a mass equal to the sum 
of the pre-encounter core of the (sub)giant and the white- 
dwarf. The newly formed white dwarf is surrounded by a 
disk formed from half the envelope of the (sub)giant before 
the encounter. If the mass of the white dwarf surpasses 
the Chandrasekhar limit, it is destroyed, without leaving 
a remnant (Nomoto & Kondo 1991 and Livio & Truran 
1985). 


3.3.4. Neutron-star or black-hole primary 

All encounters with a neutron star or black hole primary 
lead to the formation of a massive disk around the com¬ 
pact star. If the compact star had a disk prior to the col¬ 
lision, this disk is expelled. This disk accretes onto the 
compact star at a rate of 5% of the Eddington limit. An ac¬ 
creting neutron star turns into a millisecond radio-pulsar, 
or - when its mass exceeds 2M0 - into a black hole. Mu¬ 
tual encounters between neutron stars and black holes are 
ignored, as are collisions between these stars and white 
dwarfs. 


3.4- Monte Carlo Simulations of Ensembles of Encounters 

Each star in our model can encounter any of the other 
stars. To reduce computational cost, we bin the stars in 
intervals of mass and radius, and compute the probability 
for encounters between bins, giving all stars in one bin the 
same mass and radius, and, through Eq. 2, choosing their 
velocities from the same distribution. The cross section (7^ 
for a encounter with a distance of closest approach within 
d between a star from bin i and a star from bin j contains 
a geometrical and a gravitational focusing contribution: 


an = ndf 1 -I- 2G 


mi + mi 


( 10 ) 


where Vij is the relative velocity between the stars at in¬ 
finity. For the minimum separation between the two stars 
that leads to a collision d = 2{ri+ rj) is used. (The exact 
distance at which the transition between merger and bi¬ 
nary formation occurs is not known -Kochanek 1992, Lai 
et al. 1993-, we choose the factor 2 arbitrarily.) 

In the present paper, we model the stellar distributions 
as being spatially homogeneous. In order to make contact 
with astrophysical applications, we will consider our stars 
to be contained within in a fixed sphere with radius rc- 
While we can consider this radius to stand for the notion 
of ‘core radius’ in a post-collapse cluster, we want to point 
out that this interpretation is only an approximate one. 
In realistic star clusters, there is a significant drop in den¬ 
sity across the core, from the center to the core radius. 
For most stars the density drops by roughly a factor of 
three, but for the heavier stars, such as neutron stars and 
especially black holes, this factor can be much larger. 

The encounter rate Tij of stars from bin i with stars 
from bin j, anywhere in the volume of the sphere with 
radius Tc is given by two separate equations: 

r-. = /- 1) for i = j 

\ UiUj {crijVij)^Trr^ i<j, ^ 

where Ui and nj are the number densities of stars in bins i 
and j, respectively, and where () indicates averaging over 
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the the distribution of relative velocities Vij (Note that we 
should have written the last equation with an extra factor 
1/2, if we would have summed over all combinations i j, 
in order to avoid double-counting of collisions). 

Since the stars in bins i,j have Maxwellian velocity dis¬ 
tributions with root-mean-square velocity Vi and Vj , given 
by Eq.(2), the relative velocities Vij also have a Maxwellian 
distribution, with a root-mean-square velocity given by 
-I- Vj"^. Hence 




4^3 roo 
= —j= / 

v’’" Jo 


2TTd? 


l-b2G 


TTii +m 


1,2 


I^/tt 

where we have defined 
3 


r = 


2 {vi'^ + Vj'^)' 


( 12 ) 


(13) 


With this result, we write the encounter rates F^ in 
convenient units: 


„ _ rij _ Tij 

lO^pc”^ lO^pc"^ ^pc 
-3 f mi+ nij 


X [3.61 • 10' 


M, 


-b 6.31 • 10"® 

\Rq 

„-ll 


km s 


R-0 

-1 


km s 






[Myr 

The total encounter-rate follows as 

N 




ij=l 


1 


06110 


(14) 


(15) 


where N gives the total number of bins in mass and radius 
and Tenc is the average time interval between two encoun¬ 
ters. 

The stellar population in our calculation changes both 
due to encounters between stars, and due to evolution of 
the stars. The shortest evolutionary time-scale of impor¬ 
tance to us is the time scale on which the evolving stars 
expand; the fastest evolving star in the sample is used to 
set the evolution time scale 


Tev = min(i?/l?), (16) 

where R and R are the stellar radius and its time deriva¬ 
tive, respectively. 

At the beginning of each time step, we distribute the 
stars over the bins in radius and mass, calculate the num¬ 
ber densities of stars in each bin, and the evolution and 
collision time scales Tev and Tenc- The sum over all bins ij 
is less daunting as may appear at first sight, as many bins 


contain no stars. This is illustrated in Figure 3. The time 
step to be taken is then calculated as 

dt = min(0.2Tenc, Tev), (17) 

to ensure that changes in the stellar population are fol¬ 
lowed with sufficient resolution. 

At this point, a rejection technique is used to keep 
track of collisions, as follows. We choose a random number 
between 0 and 1. If this number is larger than Tdt, we 
conclude that no collision has occurred. We evolve all stars 
over a time dt, and continue with the next step. 

If the random number is smaller than Tdt, a collision 
has occurred. In calculating the sum (Eq. 15) over the 
bins, we keep track of the partial sum after addition of 
each bin combination ij. The first bin combination for 
which this growing partial sum exceeds the random num¬ 
ber identifies the bins involved in the collision. We then 
assign a sequence number to each star in bin i, and select 
one of these numbers randomly; and repeat this for bin j. 
If i and j are identical, care is taken that the same star is 
not selected twice. From the prescriptions in the previous 
section, we decide the outcome of the collision between 
the two selected stars. 

We then select another random number between 0 and 
1, to see whether a second collision has occurred. If so, we 
determine its outcome. This procedure is repeated until a 
random number larger than T dt is found, which indicates 
that no further encounter has occurred in the time step 
under consideration. 

After each time step dt a number of stars equal to 
the number of encounters that have taken place is lost 
from the stellar system; these stars have merged into sin¬ 
gle objects. For each lost star a new star is added to the 
computation, in order to guarantee a constant number of 
stars. The mass of this ‘halo guest’ is determined by the 
present-day mass-function of the cluster. 


4. A Dynamically Evolving Salpeter Mass Func¬ 
tion 

A total of two models are computed, one with a Salpeter 
type mass function which we call model S (from Salpeter), 
and one model which we call model C (from Collapsed) 
with a mass function that is affected strongly by mass 
segregation. 

In the volume of the stellar system in model S, we 
sprinkle stars according to a Salpeter mass distribution 
between O.lM0and lOOM©. The total number of stars is 
irrelevant, since we are considering these stars to be con¬ 
tained in a laboratory-type enclosure, with a thermal dis¬ 
tribution of stellar velocities. Our choice for the ‘temper¬ 
ature’ of this distribution is fixed by requiring that stars 
with a mass m = IMq will have a one-dimensional veloc¬ 
ity dispersion of 10.0 km/s, in conformity with the same 
choice made in §5. The radius of the core was chosen to be 



S.F. Portegies Zwart, P. Hut & F. Verbunt 


9 


[gs, ns 


{ms, ms\ 

{ms, gs\ 

{ms, 

wd\ 


{gs, ivd] ' 

1 

1 

1 

1 

1 

1 

1 

1 

ms 

Bss 

Yss 

wd 

+ 

disc 


0 0.2 0.4 0.6 0.8 1 


Fig. 2. Relative frequencies of various types of encounters (upper panel, the curly brackets indicate the collision) and their out¬ 
comes (lower panel), for model computation S, integrated over the duration of the calculation. Abbreviations: ms main-sequence 
star, gs giant, wd white dwarf, ns neutron star, Bss blue straggler, Yss yellow straggler. 


rc = 4 pc and the computation is started at t = 0 and ter¬ 
minated at an age of 16Gyr. For the computation of the 
encounter rate a total of 30 bins in mass, equally spaced 
in the logarithm of the mass between 0.1 and 100 M©, 
and 30 bins in radius, equally spaced in the logarithm of 
radius between 0.1 and 2000 R© are used. An additional 
bin with zero radius is used for the compact stars, i.e. the 
white dwarfs, neutron stars and black holes. 

Figure 1 shows for model computation S, the relative 
probabilities of encounters with various types of stars for 
a IM©, IR© star, at an age of the cluster of 12 Gyr. Due to 
the small encounter frequency hardly any collision prod¬ 
ucts are present in the stellar system. Only a small num¬ 
ber of blue stragglers (stars with a mass larger than the 
turnoff and with similar radii) have finite probability to 
be involved in an encounter. The most probable partners 
for an encounter with a IM©, IR© star are the stars at 
the low end of the main sequence. 

In Figure 2 we show the relative frequencies of encoun¬ 
ters of different types, and of the resulting collision prod¬ 
ucts for model S. Because the steep mass function the 
collisions rate is dominated by main-sequence stars; the 
fraction of collisions involving giants is only small. The 
most frequent type of encounter is one involving two main- 
sequence stars, leading to a main-sequence merger rem¬ 
nant with a mass smaller than the turnoff mass or a blue 
straggler when the mass of the merger exceeds the turnoff 
mass. If the mass of the merger is less than the turnoff 
mass, the product is a main-sequence star which is younger 
than primordial main-sequence stars with the same mass. 
Such a star will be left behind as a blue straggler once the 


primordial main-sequence stars leave the main-sequence. 
Yellow stragglers, i.e. giants not on the main (sub)giant 
branch of the cluster (which approximately coincides with 
the evolutionary track of a star with the turn off mass), 
can be formed directly from encounters between a main- 
sequence star and a giant, between a main-sequence star 
and a white dwarf and between a giant and a white dwarf, 
in decreasing order of importance; encounters between two 
giants are extremely rare. Our prescriptions put every 
merger product on the evolutionary track of an ordinary 
star; the presence of yellow stragglers in our calculations is 
therefore only due to the formation of giants with a mass 
larger than the turnoff mass. 

5. A More Realistic Mass Function 

The initial conditions for the mass function of the compu¬ 
tation of model C (for collapsed cluster core) are chosen 
to be more realistic, in the sense that the mass function is 
flattened due to mass segregation in the previous evolution 
of the stellar system. The lack of detailed computations 
concerning the present-day mass function in the cores of 
globular clusters, justifies our choice to use a mass function 
similar to the one described by Verbunt & Meylan (1988). 
For the mass function of model C we consider three classes 
of objects: non-degenerate stars (main-sequence stars and 
giants), white dwarfs, and neutron stars. The more mas¬ 
sive stars have all evolved, and left inert remnants (white 
dwarfs or neutron stars). We assign a certain fraction of 
the total number of stars in the stellar system to each of 
these classes. All neutron stars (5% of the total number 
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Table 2. Parameters of the different model computations, and corresponding characteristics. Subsequent columns give the name 
of the model, indication whether a Salpeter mass function or a mass function that is affected by mass segregation is used, core 
radius, 3-dimensional velocity dispersion v© time at which encounters are started, the age of the population at the start of the 
dynamical interactions, the central stellar number-density in the core, the ratio of the number of stars in the computation to 
the actual number of stars in the core, the number of encounters during the calculation per star, and the average time between 
two encounters, anywhere in the core. 


Model 

mf 

Tc 

V0 

icc 

log Uc 

/c 

^enc 

T'enc 



[pc] 

[km/s] 

[Gyr] 

[* pc-1 



[Myr] 

S 

Salpeter 

4.0 

17.3 

0 

3.93 

0.298 

0.002 

5.17 

C 

Segregated 

0.1 

17.3 

10 

6.64 

8.750 

0.660 

1.28 



Fig. 1. Relative encounter probabilities in model calculation S, 
at time t — 12Gyr, when the turnoff mass is Mto = 0.91 M©, 
for a single star with IM© and IR© as a function of mass 
and radius of the other star involved in the encounter. Darker 
shades indicate higher probabilities. The compact stars (nom¬ 
inally with zero radius) are shown as a bar below 0.1 R©: neu¬ 
tron stars between 1.34 and 2 M© and white dwarfs at lower 
masses. All other stars with radius in excess of the radius at 
the turnoff are the evolved stars. The masses of these stars is 
similar to the turnoff mass. A small fraction of blue stragglers 
is visible as an extension of the main-sequence (to the right of 
the turnoff). The vertical bar in the upper left corner presents a 
scaling to the gray shades. The lowest square corresponds to an 
encounter rate of once every 12.7 Tyr decreasing with a factor 
of two for each subsequent square. The integrated encounter 
frequency of the IM©, IR© star is 1 encounter every 1.68 Tyr. 
Almost 13% of the encounters occur with a main-sequence star 
with a mass of about ~ 0.13M© (black squares below and to 
the left). 


of stars) are assumed to have the same mass (of l.SdM©). 
The mass distribution within the two other classes are de¬ 
scribed with power-laws with a slope of a = 0 for the 
main-sequence stars and the (sub)giants and a slope of 
Of = 1 for the white dwarf progenitors. At the start of 
the dynamical modeling a total number fraction of main- 
sequence stars and giants of 70% is chosen, this number 
decreases as the stellar system evolves. The minimum ini¬ 


tial mass of a main-sequence star is chosen to be O.2M0 
instead of the O.IM© for models S. 

The numbers of stars in the different classes change as 
time evolves due to stellar evolution, encounters between 
stars, and due to the addition of a star, each time that the 
number of stars has decreased by one in a merger process. 


100 


10 

1 


0.1 


Fig. 3. Relative encounter probabilities in model calculation 
C, at time t = 12 Gyr, when the turnoff mass is Mto =0.91 M©, 
for a single star with IM© and IR© as a function of mass and 
radius of the other star involved in the encounter (similar to 
Fig. 3). The high encounter-rate and different mass function 
result in an enormous enrichment of collision products in the 
stellar system. Besides the small fraction of black holes (nomi¬ 
nally with zero radius an with a mass larger than 2M©), there 
is also a rich population of blue stragglers (in the area with a 
mass larger than the turn oft and a radius larger than about 
one R©) and yellow stragglers (stars with a radius larger than 
that of the blue stragglers). Except for the neutron stars and 
black holes (nominally with zero radius an with a mass larger 
than 2M©) all stars with mass in excess of the turnoff mass are 
the products of previous encounters. The vertical bar in the up¬ 
per left corner presents a scaling to the gray shades. The lowest 
square corresponds to an encounter rate of once every 21 Gyr 
decreasing with a factor of two for each subsequent square. 
The integrated encounter frequency of the IM©, IR© star is 
1 encounter every 3.1 Gyr. Almost 15% of the encounters oc¬ 
cur with a white dwarf with a mass of about ~ 0.7M© (black 
square in the middle and below). 
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Fig. 4. Relative frequencies of various types of encounters (upper panel) and their outcomes (lower panel), for model computation 
C, integrated over the duration of the calculation (see also Fig. 1). Apart from the variation in the relative encounter frequencies 
between various types of stars, encounters between two giants become noticeable and Thorne Zytkow objects (TZO) appear (in 
very small numbers) as the result of a collision. 


Model C has a core radius of Tc = 0.1 pc and a 1- 
dimensional velocity-dispersion for a 1 MQStar of 10 km/s. 
We switch-on the dynamics at tec =10 Gyr and terminate 
the model at t =16 Gyr. 

The number of stars used in the computation is higher 
than the calculated number of stars in the core for the 
parameters of model (7; as a result the Poissonian noise 
in our calculation is smaller than it would be in an actual 
core. 

Figure 3 shows for model computation C, the relative 
probabilities of encounters with various types of stars for 
a single IM©, IRq star, at an age of the cluster of 12 Gyr. 
At this age, products of previous encounters are already 
present in the cluster, and have a finite probability of un¬ 
dergoing another encounter. However, the most probable 
partner for an encounter with a IMq, IRq star is a white 
dwarf with a mass of about O.TMq. 

The relative importance of the various types of en¬ 
counters is very different in model C compared to model 
S, as illustrated in Figure 4, and consequently the relative 
frequencies of merger outcomes are very different as well. 
The fraction of collisions that directly result in the forma¬ 
tion of a blue straggler rises sharply as does the relative 
formation-rate of yellow stragglers and white dwarfs with 
a massive disc. Because the mass function in model C is 
flat, the region of the main sequence around the turn-off is 
well populated with massive main-sequence stars and con¬ 
sequently the total number of giants is much larger than 
in model S where a steep mass function is used. 


5.1. An Evolved H-R Diagram 



Fig. 5. Hertzsprung-Russell diagram of model C, at ca. 
t — 12 Gyr. 10^ stars (corresponding to about the total num¬ 
ber of stars in the core) were selected randomly from all stars 
involved in the simulation. 

A Hertzsprung-Russell diagram of model C after about 
12 Gyr is shown in Fig. 5. The dots (representing indi¬ 
vidual stars) that are positioned in the color magnitude 
diagram at a position that deviates from the isochrone 
of the stellar system are the result of a collision. Blue 
stragglers can be identified close to the zero-age main- 
sequence but are bluer and more luminous than the turn¬ 
off, whereas yellow stragglers are situated above the giant 
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yellow stragglers is roughly constant from the beginning 
of the dynamical simulation. 


0.02 b / 


0.0/ b / 


Fig. 6. Fraction of stars in the computation of model C that 
are blue stragglers (upper solid line) and the fraction of stars 
on the main-sequence that were left behind as blue stragglers 
when primordial stars of equal mass evolved into giants (lower 
solid line), as a function of time. Due to the slow evolution on 
the main sequence, the lower line is less susceptible to Poisso- 
nian fluctuations. The dotted lines show the fraction of stars 
that are yellow stragglers, for all yellow stragglers (upper dot¬ 
ted line) and for those that evolved from blue stragglers (lower 
dotted line). 



Fig. 7. Total number of stars formed in model C during the 
computation with dynamical encounters divided by the number 
formed from a non-dynamical model as a function of time. The 
dotted line indicates the fraction of main-sequence stars, the 
dashed line the stars on the Hertzsprung gap and (sub)giant 
branch, and the solid line the relative fraction of horizon¬ 
tal-branch stars (averaged in 100 Myr intervals). Stars on the 
(sub)giant branch are more depleted than the main-sequence 
stars as time evolves. The fraction of horizontal-branch stars 


branch. Because the stars in our calculation evolve, the 
number of collision products present at any time in the 
core is not at all proportional to their formation rate. For 
example, blue stragglers (a main sequence star with mass 
M > Mto), formed by merging of two main-sequence stars, 
often evolve into giants before our calculation is stopped, 
because of the short main-sequence lifetime of more mas¬ 
sive stars. Evolving blue stragglers turn into yellow strag¬ 
glers, and in fact most of the yellow stragglers present 
in the cluster have evolved from blue stragglers. The yel¬ 
low stragglers formed directly from collisions with giants 
evolve too fast to contribute as strongly to the presence of 
yellow stragglers. This is illustrated in Figure 6, which also 
shows that the fraction of stars that are yellow stragglers 
is rather constant throughout the computation. 

Merged main-sequence stars with a mass smaller than 
the turnoff mass upon formation are left behind as blue 
stragglers when the equally massive primordial stars leave 
the main sequence. As illustrated in Figure 6 (lower solid 
line), the fraction of such blue stragglers is relatively small. 
On the other hand, the fraction of stars that are blue 
stragglers rises rapidly at first, but levels off when the 
evolution rate of blue stragglers into yellow stragglers and 
beyond becomes competitive with their formation rate. 
Thus, the fraction of stars that are blue stragglers does 
not rise much above 3% at any given time, even though 
26% of the stars in the computation is directly turned into 
a blue straggler at some time or after after a collision. The 
dotted line in Figure 6 illustrates that the total number of 


is roughly twice as large in the stellar system where collisions 
are included. The Poissonian noise for the main-sequence stars 
is smallest, as expected, followed by that in the number of 
(sub)giants. The noise in the fraction of horizontal branch stars 
is largest. The 3 (t error bar (lower left) indicates the Poissonian 
error for the giants (the dashed line). 

Giants which undergo a collision become more mas¬ 
sive in our prescription, and thus evolve faster than their 
unperturbed counterparts. As a result, the number of gi¬ 
ants in the model is smaller than it would have been in 
a cluster without collisions, as illustrated in Figure 7. At 
the end of the computation the number of giants is de¬ 
pleted by roughly 70%. The fraction of stars on the hori¬ 
zontal branch is roughly twice that expected from a non- 
dynamically evolving stellar system. This enhancement of 
the fraction of horizontal branch stars is the result of two 
effects: most collisions between a giant and another star 
result in aging of the giant which is then evolved closer to¬ 
wards the horizontal branch and the majority of the col¬ 
lisions between a main-sequence star and a white dwarf 
results in the formation of a star that is about to termi¬ 
nate its giant lifetime (e.g. close to or on the horizontal 
branch). 

6. Conclusions 

The models discussed in this paper are very crude in their 
treatment of the encounter processes, of the result of a 
collision between two stars, and of the evolution of the 
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merger products. Apart from these approximations and 
the fact that we use a stellar evolution model for popula¬ 
tion I instead of pop. II stars, the adopted mass function 
is also highly uncertain. Nonetheless, some interesting re¬ 
sults can be delineated. 

Comparison with the calculations of Davies & Benz 
(1995) shows the effect of allowing the merger products 
to evolve. An immediate consequence of this is the lower 
prediction for the number of blue and yellow stragglers 
present in the cluster (as is clear from Figure 6). The for¬ 
mation rates of blue and yellow stragglers give a poor in¬ 
dication for the actual number of stragglers present in the 
cluster at a particular instant. 

Due to the low density of model S the collision fre¬ 
quency is small. The steep Salpeter mass function also 
suppresses the encounter rate and the production of stel¬ 
lar curiosities; the majority of the collisions involve two 
rather low mass main-sequence stars which results in a 
merger that evolves too slow to produce a blue stragglers 
within the time span of the simulation. 

The Hertzsprung-Russell diagram of our model cluster 
(model C) shows that blue stragglers close to the turn-off 
point lie on the main sequence, whereas blue stragglers 
above the turnoff point are mostly found at some distance 
from the main sequence. The reason for this is that colli¬ 
sions only become important in the cluster when an ini¬ 
tial period of low density is followed by the contraction 
of the cluster core. The more massive blue stragglers are 
formed in collisions between stars close to the terminal- 
age main sequence, and evolve relatively quickly. Blue 
stragglers close to the turn-off are formed in collisions be¬ 
tween relatively low-mass stars which did not evolve very 
far away from the zero-age main sequence, and therefore 
also the merger products are close to the zero-age main 
sequence, and evolve slowly. Thus, the point where blue 
stragglers have left the main sequence gives an indication 
of the time when collisions in the cluster became frequent 
(see also Portegies Zwart 996a). 

Our model C predicts a depletion of giants, in the core 
only, up to ~ 50% shortly after tec relative to a collision¬ 
less stellar system, in globular clusters with a collapsed 
core where the fraction of horizontal-branch stars is en¬ 
hanced. Consequently the depletion of giants relative to 
the number of horizontal branch stars is strongly present 
in the high-density stellar system. Collisions between sin¬ 
gle stars cannot explain the observation that giants can be 
depleted well outside the core or completely absent in it, 
as observed in the core of M 15 (Djorgovski et al. 1991). 

In our simulated cluster cores the total number of 
white dwarfs that exceed the Chandrasekhar limit due to 
accretion from a circum-stellar disc is small, even in the 
cluster simulation with the highest density. In model C 
8% of the white dwarfs experience an accretion-induced 
collapse, which (after correction for the ratio /c between 
the number of stars in the model and in an actual core 
- see Table 2) corresponds to 190 supernovae of type la 


during the 6 Gyr of our calculation. If all of these collapses 
would lead to the formation of a neutron star, and if all of 
these would remain in the core, this would be a substan¬ 
tial addition to the total number of neutron stars in the 
core, which is about 460 (after correction for fc) at the 
start of our calculation. This result, however, strongly de¬ 
pends on the adopted mass function for the white dwarfs. 
The formation-rate of neutron stars with an accretion disc 
and the subsequent formation of a recycled pulsar or black 
hole is (to first order) linearly dependent on the number of 
neutron stars, which depends not only on the initial mass 
function but also on the subsequent mass segregation in 
the cluster. 

The encounter rates between neutron stars and main- 
sequence stars are similar in our calculations to the rates 
found in the calculations by Verbunt & Meylan (1988), 
by Di Stefano & Rappaport (1992) and by Davies & Benz 
(1995). After a collision between a neutron star and an¬ 
other cluster member the merged object becomes visible 
as an X-ray source (for at most 1 Gyr) after which it be¬ 
comes a recycled pulsar or, if its mass exceeds 2 Mq, a 
black hole. The total number of such X-ray sources, re¬ 
cycled pulsars or black holes scales linearly, in first order, 
with the number of neutron stars in the cluster core, which 
is rather uncertain. 

Our computations reveal that collisions between single 
stars result in a small number of recycled pulsars: about 70 
are formed in a core (after correction for fc) according to 
model C. Whether this is enough to explain the observed 
numbers is not clear. The intrinsic luminosity distribution 
of millisecond pulsars, and hence the fraction of them that 
is detectable in a typical cluster, is not known; and some 
clusters with high encounter rates show remarkably few 
recycled pulsars, the globular cluster NGC 6342 is an ex¬ 
ample (see Lyne 1993). 
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